

      program main

      use globais
      implicit none

      integer :: simula
      integer :: ierr
      integer :: it!< contador de instantes de tempo
      integer :: uerros =26!< unidade para impressão dos erros do modelo em relação ao modelo real 

      double precision, dimension(:), allocatable :: Staux !<variavel auxiliar para verificacao da transformada inversa de fourier
    
      character*5 :: filesism
      integer :: iprint!< indica se na simulação as chamadas ao pswigp serão plotadas na tela
      double precision :: freqmax !< frequência máxima 


      call entrada(simula)

!       freqmax=20.0d0
!       call calcMisfitFinal(gruposSR,Nsource,Nrec,Nttotal,ntps,uerros,freqmax)! cálculo do residuo no tempo e na frequência, a partir dos sismogramas salvos
!       stop

!       open(1,file='outfonte.txt')
!       write(1,*)'#tempo, parte real, parte imaginaria '
!       Nttotalaux = 2*Nomega+1
!       dtaux = tf/dfloat(Nttotalaux-1)
!       do it=0,Nttotalaux-1
! 	  Staux = 0.
! 	  do kw=1,2*Nomega
! 	      Staux = Staux + Sw(kw)/tf * EXP(imag*omegas(kw)*it*dtaux)
! 	  enddo
! 	  write(1,*)it*dtaux, real(Staux),aimag(Staux)
!       enddo
!       close(1)
	



!       call fonte() !< transformada de Fourier da fonte


      if(simula.eq.0)then

	  allocate(St(nttotal),stat=ierr)
	  call rickerwavelet(tf,nttotal,fcorte,fpico,ampS,tp,St)

	  allocate(Sw(Nomega),stat=ierr)!< transformada de Fourier
	  call dft(Sw,omegas,Nomega,ntps,St,nttotal)

	  open(1,file='outfontefreq.txt')
	  do kw=1,nomega
	     write(1,*)omegas(kw)/(2.0d0*Pi), real(Sw(kw)),aimag(Sw(kw)),sqrt(real(Sw(kw))**2+aimag(Sw(kw))**2) 
	  enddo
	  close(1)

	  call inversao()

	  deallocate(Sw)


	  open(uerros,file='erromodel.txt')

	  call calcErro(nx,nz,npml,uerros)


	  filesism = 'sismA' !< sismogramas do modelo Aproximado
	  iprint = 0 !< não plota as chamadas ao pswigp na tela
	  call simulacao(2,filesism,iprint)!simulação no domínio do tempo, salvando os sismogramas

	  freqmax = omegas(nomega)/(2.0d0*Pi)
	  call calcMisfitFinal(gruposSR,Nsource,Nrec,Nttotal,ntps,uerros,freqmax)! cálculo do residuo no tempo e na frequência, a partir dos sismogramas salvos



	  close(uerros)	  

      elseif(simula.gt.0)then

	  if(simula.eq.1)then !<domínio da frequência

	      allocate(St(nttotal),stat=ierr)
	      call rickerwavelet(tf,nttotal,fcorte,fpico,ampS,tp,St)

	      allocate(Sw(2*Nomega),stat=ierr)!< transformada de fourier
	      call dft(Sw,omegas,Nomega,ntps,St,nttotal)
	      do kw=1,Nomega
		  Sw(2*Nomega-(kw-1))=CONJG(Sw(kw))
	      enddo
	      open(1,file='outfontefreq.txt')
	      write(1,*)'#freq em Hz, parte real, parte imaginaria, amplitude '
	      do kw=1,2*Nomega
		  write(1,*)omegas(kw)/(2.*Pi), real(Sw(kw)),aimag(Sw(kw)),sqrt(real(Sw(kw))**2+aimag(Sw(kw))**2)
	      enddo
	      close(1)


	      allocate(Staux(nttotal),stat=ierr)!< verificacao da transformada de fourier via transfomada inversa
	      call invdft(Sw,omegas,Nomega,ntps,Staux,nttotal)
	      open(1,file='outfonte.txt')
	      write(1,*)'#tempo, parte real '
	      do it=0,Nttotal-1
		  write(1,*)it*tf/dfloat(nttotal),Staux(it+1)
	      enddo
	      close(1)
	      deallocate(Staux)!<fim da verificacao

	  elseif(simula.eq.2)then

	      !vou gerar a fonte só depois que verificar deltat para estabilidade
! 	      allocate(St(nttotal),stat=ierr)
! 	      call rickerwavelet(tf,nttotal,fcorte,fpico,ampS,tp,St)

	  endif


	  call modelo2dcamadas()

	  filesism = 'sismD' !< geração dos sismogramas que serão usados como Dados para a inversão 
	  iprint = 1
	  call simulacao(simula,filesism,iprint)!simulação no domínio da frequencia (simula=1) ou no domínio do tempo (simula=2)

	  if(simula.eq.1)then !<domínio da frequência
	      deallocate(Sw)
	  endif


      endif 






      deallocate(csix)
      deallocate(csiz)
      deallocate(c)
      deallocate(creal4)
      deallocate(rho)
      deallocate(b)
      deallocate(K)
      deallocate(Preceiv)

      deallocate(St)
      deallocate(omegas)
      deallocate(ireceiv)
      deallocate(jreceiv)


      stop


      end program main